# Loughney et al., "Tectonic influence on Cenozoic mammal richness and sedimentation history of the Basin and Range, western North America"
# Science Advances
# Partial-correlation analyses and output

# load required package
library(ppcor)

# load files
macrostrat <- read.csv("macrostrat_June-2021.csv", header = TRUE, stringsAsFactors = FALSE)
fossils <- read.csv("fossiliferous_units_January-2021.csv", header = TRUE, stringsAsFactors = FALSE)
deformation <- read.csv("deformation_rates_January-2021.csv", header = TRUE, stringsAsFactors = FALSE)
area <- read.csv("area_change.csv", header = TRUE, stringsAsFactors = FALSE)

# calculate SAR of fossiliferous units of the Basin and Range
> brfossilsar <- fossils$NB_SAR + fossils$CB_SAR + fossils$SB_SAR

# partial-correlation analyses
# ------------------------------------------------------------------------------------------------------------------------
# species richness - deformation rates (area-change rates controlled) to 0 Ma
> pcor.test(brAllRichness[,2], deformation$BR_rate, area$BR_rate_change_sum, method = "spearman")
     estimate   p.value  statistic  n gp   Method
1 -0.01916207 0.8739754 -0.1592013 72  1 spearman

# species richness - deformation rates (area-change rates controlled) to 0.5 Ma
> pcor.test(brAllRichness[2:72,2], deformation$BR_rate[2:72], area$BR_rate_change_sum[2:72], method = "spearman")
    estimate   p.value statistic  n gp   Method
1 0.04115573 0.7351545 0.3396667 71  1 spearman

# species richness - area-change rates (deformation rates controlled) to 0 Ma
> pcor.test(brAllRichness[,2], area$BR_rate_change_sum, deformation$BR_rate, method = "spearman")
  estimate      p.value statistic  n gp   Method
1 0.670074 1.649519e-10  7.498444 72  1 spearman

# species richness - area-change rates (deformation rates controlled) to 0.5 Ma
> pcor.test(brAllRichness[2:72,2], area$BR_rate_change_sum[2:72], deformation$BR_rate[2:72], method = "spearman")
   estimate      p.value statistic  n gp   Method
1 0.6415787 2.170995e-09  6.897261 71  1 spearman

# macrostrat SAR - deformation rates (area-change rates controlled) to 0 Ma
> pcor.test(macrostrat$BR_SAR, deformation$BR_rate, area$BR_rate_change_sum, method = "spearman")
    estimate    p.value statistic  n gp   Method
1 -0.2741787 0.02067896 -2.368254 72  1 spearman

# macrostrat SAR - area-change rates (deformation rates controlled) to 0 Ma
> pcor.test(macrostrat$BR_SAR, area$BR_rate_change_sum, deformation$BR_rate, method = "spearman")
   estimate      p.value statistic  n gp   Method
1 0.6895057 2.947489e-11  7.907801 72  1 spearman

# macrostrat package thickness - deformation rates (area-change rates controlled) to 0 Ma
> pcor.test(macrostrat$BR_thickness, deformation$BR_rate, area$BR_rate_change_sum, method = "spearman")
    estimate   p.value statistic  n gp   Method
1 -0.1659613 0.1666015 -1.397965 72  1 spearman

# macrostrat package thickness - area-change rates (deformation rates controlled) to 0 Ma
> pcor.test(macrostrat$BR_thickness, area$BR_rate_change_sum, deformation$BR_rate, method = "spearman")
   estimate      p.value statistic  n gp   Method
1 0.6687137 1.851949e-10  7.470886 72  1 spearman

# fossiliferous unit SAR - deformation rates (area-change rates controlled) to 0 Ma
> pcor.test(brfossilsar, deformation$BR_rate, area$BR_rate_change_sum, method = "spearman")
   estimate      p.value statistic  n gp   Method
1 0.4190172 0.0002759193  3.833372 72  1 spearman

# fossiliferous unit SAR - area-change change rates (deformation rates controlled) to 0 Ma
> pcor.test(brfossilsar, area$BR_rate_change_sum, deformation$BR_rate, method = "spearman")
   estimate      p.value statistic  n gp   Method
1 0.3922419 0.0007167579  3.542059 72  1 spearman

# thickness of fossiliferous units - deformation rates (area-change rates controlled) to 0 Ma
> pcor.test(fossils$BR_thickness, deformation$BR_rate, area$BR_rate_change_sum, method = "spearman")
   estimate      p.value statistic  n gp   Method
1 0.4859873 1.735697e-05  4.619074 72  1 spearman

# thickness of fossiliferous units - area-change rates (deformation rates controlled) to 0 Ma
> pcor.test(fossils$BR_thickness, area$BR_rate_change_sum, deformation$BR_rate, method = "spearman")
   estimate     p.value statistic  n gp   Method
1 0.3157603 0.007309222  2.764327 72  1 spearman
